function pairPulseAnalysis(handles)
%  Analyze pair pulse map data.
%   pairPulseAnalysis(handles)
%
%   See also

% Author: gs
% Created: Dec 2006
% Copyright 2006
% -------------------------------------
if nargin == 0
    try
        handles = evalin('base', 'handles');
    catch
        disp('Unable to find variable ''handles'' in workspace');
    end
end

% FIGURE
hfig = figure('Color', 'w', 'DoubleBuffer', 'on', 'Units', 'normalized', ...
    'Position', [0.1644    0.2900    0.6756    0.5000], ...
    'Name', 'pairPulseAnalysis');
if strcmp(handles.data.analysis.eventPolaritySyn, 'down')
    plotPolarityFactor = 1;
    colormap(flipud(colormap(jet2(256))));
    eventThresh = -2; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <<<<<<< set this
elseif strcmp(handles.data.analysis.eventPolaritySyn, 'up')
    plotPolarityFactor = -1;
    colormap(jet2(256));
    eventThresh = 1.5; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <<<<<<< set this
end

% settings
PPRclim = [0 2]; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <<<<<<< set this
maskFlag = 1;

% ANALYSES ================================
% mean
M1 = handles.data.map.mapStack.mapAvgMean;
M2 = handles.data.map.mapStack.mapAvgMean4th;
PPRmean = M2./M1;
if plotPolarityFactor == 1 % downgoing events
    Mask = double(M1>eventThresh);
else % upgoing events
    Mask = double(M1>eventThresh);
end
Mask(Mask==0) = NaN;
M1meanMask = Mask;
% mean
M1 = handles.data.map.mapStack.mapAvgMin;
M2 = handles.data.map.mapStack.mapAvgMin4th;
PPRmin = M2./M1;
if plotPolarityFactor == 1 % downgoing events
    Mask = double(M1>eventThresh);
else % upgoing events
    Mask = double(M1>eventThresh);
end
Mask(Mask==0) = NaN;
M1minMask = Mask;
% mean
M1 = handles.data.map.mapStack.mapAvgMax;
M2 = handles.data.map.mapStack.mapAvgMax4th;
PPRmax = M2./M1;
if plotPolarityFactor == 1 % downgoing events
    Mask = double(M1>eventThresh);
else % upgoing events
    Mask = double(M1>eventThresh);
end
Mask(Mask==0) = NaN;
M1maxMask = Mask;

% ==========================================

% MEAN maps
map = handles.data.map.mapStack.mapAvgMean;
if maskFlag
    map = map.*M1meanMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub(1) = subplot(3,3,1);
imagesc(map);
daspect([1 1 1]);
title('Pulse#1: mean');
clims = get(gca, 'CLim');
%
map = handles.data.map.mapStack.mapAvgMean4th;
if maskFlag
    map = map.*M1meanMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub(2) = subplot(3,3,2);
imagesc(map);
daspect([1 1 1]);
clims = [clims get(gca, 'CLim')];
set(hsub(1:2), 'CLim', [min(clims) max(clims)]);
colorbar;
title('Pulse#2: mean');
%
map = PPRmean;
if maskFlag
    map = map.*M1meanMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub = subplot(3,3,3);
imagesc(map);
daspect([1 1 1]);
set(gca, 'CLim', [PPRclim(1) PPRclim(2)]);
colorbar;
title('PPR: mean');

% MIN maps
map = handles.data.map.mapStack.mapAvgMin;
if maskFlag
    map = map.*M1minMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub(3) = subplot(3,3,4);
imagesc(map);
daspect([1 1 1]);
clims = get(gca, 'CLim');
title('Pulse#1: min');
%
map = handles.data.map.mapStack.mapAvgMin4th;
if maskFlag
    map = map.*M1minMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub(4) = subplot(3,3,5);
imagesc(map);
daspect([1 1 1]);
clims = [clims get(gca, 'CLim')];
set(hsub(3:4), 'CLim', [min(clims) max(clims)]);
colorbar;
title('Pulse#2: min');
%
map = PPRmin;
if maskFlag
    map = map.*M1minMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub = subplot(3,3,6);
imagesc(map);
daspect([1 1 1]);
set(gca, 'CLim', [PPRclim(1) PPRclim(2)]);
colorbar;
title('PPR: min');

% MAX maps
map = handles.data.map.mapStack.mapAvgMax;
if maskFlag
    map = map.*M1maxMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub(5) = subplot(3,3,7);
imagesc(map);
daspect([1 1 1]);
clims = get(gca, 'CLim');
title('Pulse#1: max');
%
map = handles.data.map.mapStack.mapAvgMax4th;
if maskFlag
    map = map.*M1maxMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub(6) = subplot(3,3,8);
imagesc(map);
daspect([1 1 1]);
clims = [clims get(gca, 'CLim')];
set(hsub(5:6), 'CLim', [min(clims) max(clims)]);
colorbar;
title('Pulse#2: max');
%
map = PPRmax;
if maskFlag
    map = map.*M1maxMask;
end
map(isnan(map)) = plotPolarityFactor * inf;
hsub = subplot(3,3,9);
imagesc(map);
daspect([1 1 1]);
set(gca, 'CLim', [PPRclim(1) PPRclim(2)]);
colorbar;
title('PPR: max');
